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Abstract 

1/1/e are interested in the strong convergence and almost sure stability of Euler-Maruyama (EM) type 
approximations to the solutions of stochastic differential equations (SDEs) with non-linear and non- 
Lipschitzian coefficients. Motivation comes from finance and biology where many widely applied models 
do not satisfy the standard assumptions required for the strong convergence. In addition we examine the 
globally almost surely asymptotic stability in this non-linear setting for EM type schemes. In particular, we 
present a stochastic counterpart of the discrete LaSalle principle from which we deduce stability properties 
for numerical methods. 

Key words: Dissipative model, super-linear growth, stochastic differential equation, strong conver- 
gence, backward Euler-Maruyama scheme, implicit method, LaSalle principle, non-linear stability, almost 
sure stability. 

AMS Subject Clasification: 65C30, 65L20, 60H10 

1 Introduction 

Throughout this paper, let {VL,T ^{Tt\t>o,^) be a complete probability space with a filtration {Tt\t>a 
satisfying the usual conditions (that is to say, it is right continuous and increasing while J^o contains all P- 
nuU sets). Let w{t) — {wi{t), Wd{t))'^ be a c?-dimensional Brownian motion defined on the probability 
space, where T denotes the transpose of a vector or a matrix. In this paper we study the numerical 
approximation of the stochastic differential equations (SDEs) 

dx{t) = f{x{t))dt + g{x{t))dw{t). (1.1) 

Here x{t) e R" for each t > and / : M" M" and 5 : M" ^- R"^''. For simplicity we assume that 
xq G M" is a deterministic vector. Although the method of Lyapunov functions allows us to show that 
there are solutions to a very wide family of SDEs (see e.g. [MJIM]), in general, both the explicit solutions 
and the probability distributions of the solutions are not known. We therefore consider computable 

'Department of Mathematics and Statistics, University of Strathclyde, Glasgow, Gl IXH, Scotland, UK 
(x.mao@strath.ac.uk). 

t Mathematical Institute, University of Oxford, Oxford, OXl 3LB, UK (lukas.szpruchamaths.ox.ac.uk). 



discrete approximations that, for example, could be used in Monte Carlo simulations. Convergence and 
stability of these methods are well understood for SDEs with Lipschitz continuous coefficients: see [26] 
for example. Our primary objective is to study classical strong convergence and stability questions for 
numerical approximations in the case where / and g are not globally Lipschitz continuous. A good 
motivation for our work is an instructive conditional result of Higham et al. |18j . Under the local 
Lipschitz condition, they proved that uniform boundedness of moments of both the solution to 
and its approximation are sufficient for strong convergence. That immediately raises the question of 
what type of conditions can guarantee such a uniform boundedness of moments. It is well known that 
the classical linear growth condition is sufficient to bound the moments for both the true solutions and 
their EM approximation [26) 134) . It is also known that when we try to bound the moment of the true 
solutions, a useful way to relax the linear growth condition is to apply the Lyapunov-function technique, 
with V{x) = ||a;|p. This leads us to the monotone condition [33]. More precisely, if there exist constants 
a, /3 > such that the coefficients of equation satisfy 

{x,f{x)) + hg{x)\\^ <a + (3\\xW' for all x G M", (1.2) 



then 



2' 



sup E||a;(t)||^ < oo VT > 0. (1.3) 

0<t<T 



Here, and throughout, ||x|| denotes both the Euclidean vector norm and the Frobenius matrix norm and 
{x, y) denotes the scalar product of vectors x,y G R". However, to the best of our knowledge there is no 
result on the moment bound for the numerical solutions of SDEs under the monotone condition p. 21) . 
Additionally, Hutzenthaler et al. [l^ proved that in the case of super-linearly growing coefficients, the 
EM approximation may not converge in the strong L^-sense nor in the weak sense to the exact solution. 
For example, let us consider a non-linear SDE 

dx{t) ^ {fi- ax{tf)dt + I3x{tfdw{t), (1.4) 

where ^i,a,P> and a > ■^f3'^. In order to approximate SDE (|1.4p numerically, for any At, we define 
the partition pAt := {tk = kAt : k = 0, 1, 2, . . . , iV} of the time interval [0,T], where N At = T and 
T > 0. Then we define the EM approximation Yt^^ ~ x{tk) of (|1.4p by 

= Yt. + - aYl)At + PY^Awt,, (1.5) 

where Awt^. = w{tk+i) — w{tk)- It was shown in [23 that 

lim E||rt„||2 = oo. 

At^O 

On the other hand, the coefficients of (|1.4p satisfy the monotone condition (|1.2p so (|1.3p holds. Hence, 
Hutzenthaler et al. 23 concluded that 

^m¥.\\x{T)-YtA^^^- 

It is now clear that to prove the strong convergence theorem under condition ()1.2p it is necessary to 
modify the EM scheme. Motivated by the existing works fTS^ and f2P we consider implicit schemes. 
These authors have demonstrated that a backward Euler-Maruyama (BEM) method strongly converges 
to the solution of the SDE with one-sided Lipschitz drift and linearly growing diffusion coefficients. So 
far, to the best of our knowledge, most of the existing results on the strong convergence for numerical 
schemes only cover SDEs where the diffusion coefficients have at most linear growth [S] ISHl HZl [221 [21] ■ 
However, the problem remains essentially unsolved for the important class of SDEs with super-linearly 
growing diffusion coefficients. We are interested in relaxing the conditions for the diffusion coefficients 
to justify Monte Carlo simulations for highly non-linear systems that arise in financial mathematics, 
[3 [HI [1 il [H [2S] , for example 

dx{t) ^ {fi- ax''{t))dt + l3xP{t)dw{t), r,p>l, (1.6) 
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where fJ,,a,[3 > 0, or in stochastic population dynamics [351 151 1551 BTl 111) , for example 

dx{t) = diag{xi{t),X2{t), x„(t))[(fe + Ax'^{t))dt + x{t)dw{t)], (1.7) 

where h = . . . , &n)"^, 2;^(i) = {x\{t), . . . , x^(i)„)-^ and matrix A = [ay]i<i j<„ is such that Xmax{A + 
A^) < 0, where Xmax{A) = sup^g^n x'^Ax. The only results we know, where the strong convergence 
of the numerical approximations was considered for super-linear diffusion, is Szpruch et al. |43| and Mao 
and Szpruch [37]. In [33] authors have considered the BEM approximation for the following scalar SDE 
which arises in finance [2], 

dx{t) — {a-ix{t)~^ — an + aix(t) — a2x{tY)dt + ax{tYdw{t) r, p > 1. 

In [37] , this analysis was extended to the multidimensional case under specific conditions for the drift and 
diffusion coefficients. In the present paper, we aim to prove strong convergence under general monotone 
condition (jl.2p in a multidimensional setting. We believe that this condition is optimal for boundedness 
of moments of the implicit schemes. The reasons that we are interested in the strong convergence are: 
a) the efficient variance reduction techniques, for example, the multi-level Monte Carlo simulations [12) . 
rely on the strong convergence properties; b) both weak convergence (261 and pathwise convergence [251 
follow automatically. 

Having established the strong convergence result we will proceed to the stability analysis of the under- 
lying numerical scheme for the non-linear SDEs (II. ip under the monotone-type condition. The main 
problem concerns the propagation of an error during the simulation of an approximate path. If the 
numerical scheme is not stable, then the simulated path may diverge substantially from the exact solu- 
tion in practical simulations. Similarly, the expectation of the functional estimated by a Monte Carlo 
simulation may be significantly different from that of the expected functional of the underlying SDEs 
due to numerical instability. Our aim here is to investigate almost surely asymptotic properties of the 
numerical schemes for SDEs (|l.ip via a stochastic version of the LaSalle principle. LaSalle, [27], improved 
significantly the Lyapunov stability method for ordinary differential equations. Namely, he developed 
methods for locating limit sets of nonautonomous systems |13l I27j . The first stochastic counterpart of 
his great achievement was established by Mao [33] under the local Lipschitz and linear growth condi- 
tions. Recently, this stochastic version was generalized by Shen et al. [42] to cover stochastic functional 
differential equations with locally Lipschitz continuous coefficients. Furthermore, it is well known that 
there exists a counterpart of the invariance principle for discrete dynamical systems [28) . However, there 
is no discrete stochastic counterpart of Mao's version of the LaSalle theorem. In this work we investigate 
a special case, with the Lyapunov function V{x) — of the LaSalle theorem. We shall show that 

the almost sure global stability can be easily deduced from our results. The primary objectives in our 
stability analysis are: 

• Ability to cover highly nonlinear SDEs; 

• Mild assumption on the time step - A(a)-stability concept [15] . 

Results which investigate stability analysis for numerical methods can be found in Higham [TB] [15] for 
the scalar linear case. Baker et al. [4] for the global Lipschitz and Higham et al. [T9| for one-sided 
Lipschitz and the linear growth condition. 

At this point, it is worth mentioning how our work compares with that of Higham et al. [18^. Theorem 
3.3 in their paper is a very important contribution to the numerical SDE theory. The authors proved 
strong convergence results for one-sided Lipschitz and the linear growth condition on drift and diffusion 
coefficients, respectively. What differentiates our work from [18] are: a) We significantly relax the 
linear growth constraint on the diffusion coefficient and we only ask for very general monotone type 
growth; b) Our analysis is based on a more widely applied BEM scheme in contrast to the split-step 
scheme introduced in their paper. An interesting alternative to the implicit schemes for numerical 
approximations of SDEs with non-globally Lipschitz drift coefficient recently appeared in 1221 . However 
the stability properties of this method are not analysed. 

In what follows, for notational simplicity, we use the convention that C represent a generic positive 
constant independent of At, the value of which may be different for different appearances. 
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The rest of the paper is arranged as follows. In section 2, we introduce the monotone condition under 
which we prove the existence of a unique solution to equation along with appropriate bounds that 

will be needed in further analysis. In Section 3 we propose the 6'-EM scheme, which is known as the BEM 
when 6* = 1, to approximate the solution of equation (jl.ip . We show that the 2nd moment of the 0-EM, 
can be bounded under the monotone condition plus some mild assumptions on / and g. In Section 4 we 
introduce a new numerical method, which we call the forward- backward Euler-Maruyama (FBEM) . The 
FBEM scheme enables us to overcome some measurability difficulties and avoid using Malliavin calculus. 
We demonstrate that both the FBEM and the d-FiM do not differ much in the L^-sense. Then we prove 
a strong convergence theorem on a compact domain that is later extended to the whole domain. We also 
perform a numerical experiment that confirms our theoretical results. Section 5 contains the stability 
analysis, where we prove a special case of the stochastic LaSalle theorem for discrete time processes. 



2 Existence and Uniqueness of Solution 

We require the coefficients / and g in (jl.ll) to be locally Lipschitz continuous and to satisfy the monotone 
condition, that is 

Assumption 2.1. Both coefficients f and g in (11.11) satisfy the following conditions: 

Local Lipschitz condition. For each integer m > 1, there is a positive constant C{m) such that 

\\f{^) - f{v)\\ + M^) - 9{V)\\ < C{m)\\x - y\\ 
for those x,y €z M" with \\x\\ V ||y|| < m. 

Monotone condition . There exist constants a and /3 such that 

{x,fix))+^\\gix)\\'<a + (3\\x\\' (2.1) 

for allxe M". 

It is a classical result that under Assumption 12. 1[ there exists a unique solution to (|l.ip for any given 
initial value a;(0) = xq G K", [inilSl]- The reason why we present the following theorem with a proof here 
is that it reveals the upper bound for the probability that the process x{t) stays on a compact domain 
for finite time T > 0. The bound will be used to derive the main convergence theorem of this paper. 

Theorem 2.2. Let Assumvtion \2. 1\ hold. Then for any given initial value x{0) ^ xq Cz K", there exists 
a unique, global solution {x{t)}t>o to equation il.l]} . Moreover, the solution has the properties that for 
any T > 0, 

E\\x{T)f < (|la;o|P + 2ar)exp(2/3T), (2.2) 

and 

P(.„<T)<«±^^^^)^^Pi«, (2.3) 
where m is any positive integer and 

r„, = inf{i>0: |la;(t)|l > m}. (2.4) 

Proof. It is well known that under Assumption 12. 1[ for any given initial value xq € M" there exists a 
unique solution x{t) to the SDE p.ip . [?, 34]. Therefore we only need to prove that (|2.2p and (|2.3I) hold. 
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Applying the Ito formula to the function V{x,t) = Ha^lP, we compute the diffusion operator 

LVix,t) = 2[{xJix)) + ^\\gix)f). 

By Assumption 12.11 

LV{x,t) <2a + 2l3\\xf. (2.5) 

Therefore ^ 

E||a;(t Ar™)||^ < ||xo||^ + 2ar+ / 2^E\\x{s A Trn)fds. 

Jo 

The Gronwall inequality gives 

E|la;(r A t„,W < (IkolP + 2aT) exp(2/3T). (2.6) 

Hence 

P(t™ < T)m2 < [\\xo\\^ + 2aT] exp(2/3T). 
Next, letting m oo in (|2.6p and applying Fatou's lemma, we obtain 

E\\x{T)f < [||xof + 2ar]exp(2/3T), 

which gives the other assertion (|2.2p and completes the proof. □ 



3 The ^-Euler-Maruyama Scheme 

As indicated in the introduction, in order to approximate the solution of we will use the d-EM 

scheme. Given any step size At, we define a partition V^t '■— {tk — kAt : fc = 0, 1, 2, ...} of the half line 
[0, oo), and define 

Xt,^, = Xt, + ef{Xt,^,)At + (1 - e)f{Xt,)At + g{XtJAwt,, (3.1) 

where Awt^. = Wt^^^ — Wt^, and Xt^ = xq. The additional parameter 9 e [0, 1] allows us to control the 
implicitness of the numerical scheme, that may lead to various asymptotic behaviours of equation p.ip . 
For technical reasons we always require 9 > 0.5. 

Since we are dealing with an implicit scheme we need to make sure that equation (13.11) has a unique 
solution Xt^^-^ given Xt,^. To prove this, in addition to Assumption 12.11 wc ask that function / satisfies 
the one-sided Lipschitz condition. 

Assumption 3.1. One-sided Lipschitz condition. There exists a constant L > such that 

{x - y, fix) - fiy)) <L\\x^ yf Vx, y G M". 

It follows from the fixed point theorem that a unique solution Xt^^-^ to equation (13.11) exists given Xt^^, 
provided At < (see [57] for more details). From now on we always assume that At < In order 
to implement numerical scheme (j3.1l) we define a function F : R" R" as 

F{x) ^x-9f{x)At. (3.2) 

Due to Assumption 13. 11 there exists an inverse function and the solution to p. II) can be represented 
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in the following form 



Xt,^, = F-\Xt, + (1 - e)fiXt,)At + giXt,)Awt,). 

Clearly, Xt^. is J^t,, -measurable. In many applications, the drift coefBcient of the SDEs has a cubic or 
quadratic form, whence the inverse function can be found explicitly. For more complicated SDEs we can 
find the inverse function using root-finding algorithms, such as Newton's method. 



3.1 Moment Properties of 6'-EM 

In this section we show that the second moment of the O-FiM (j3.1D is bounded fTheorem l3.6p . To achieve 
the bound we employ the stopping time technique, in a similar way as in the proof of Theorem 12.21 
However, in discrete time approximations for a stochastic process, the problem of overshooting the level 
where we would like to stop the process appears, [71 [51 132). 

Due to the implicitness of scheme p.ip . an additional but mild restriction on the time step appears. 
That is, from now on, we require At < At* , where At* e (0, {max{L, 2j3}9)^^) with f3 and L defined in 
Assumptions 12.11 and 13.11 respectively. 

The following lemma shows that in order to guarantee the boundedness of moments for Xt^. defined by 
(|3.ip it is enough to bound the moments of F{Xt^), where F is defined by 



Lemma 3.2. Let Assumvtion \2. 1\ hold. Then for F(x) — x ~ 6f{x)At we have 

P < (1 - 2/39At)-^ [\\Fix)f + 2eaAt] Vx G R". 



Proof. Writing = {F{x),F{x)) and using Assumption 12. II we arrive at 

\\F{x)r =\\xr - "^oix, fix)) At + e^f{x)rAt^ 

>(1 - 2(30At)\\xf - 20aAt, 
and the assertion follows. □ 

We define the stopping time Am by 

X„^^mi{k:\\XtJ\>m}. (3.3) 

We observe that when k G [0, Am(w)], ||A^tfc_i(w)|| < m, but we may have ||Xt^(a;)|| > m, so the following 
lemma is not trivial. 

Lemma 3.3. Let Assumvtions lKTi \3.1\ hold, and > 0,5. Then for p > 2 and sufficiently large integer 
m, there exists a constant C{p,m), such that 

E[||XtJ|fl[o,A„](fc)] <C(p,m) for any fc > 0. 

Proof. The proof is given in the Appendix. □ 

For completeness of the exposition we recall the discrete Gronwall inequality, that we will use in the 
proof of Theorem [ 
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Lemma 3.4 (The Discrete Gronwall Inequality). Let M he a positive integer. Let Uk and Vk be non- 
negative numbers for k=0,l,...,M. If 



fc-i 

Uk < + '^^VjUj, Vfc = 1, 2, M, 



then 

/k-i 

Uk < uo exp ^ ) , V/c = 1, 2, M. 
\j=o 

The proof can be found in Mao et al. [35]. To prove the boundedness of the second moment for the 
^-EM (j3.ip . we need an additional but mild assumption on the coefficients / and g. 

Assumption 3.5. The coefficients of equation lil.l]) satisfy the polynomial growth condition. That is, 
there exists a pair of constants h > I and C{h) > such that 

\\f{x)\\y\\g{x)\\<C{h){l + \\xt), VxeR". (3.4) 



Let us begin to establish the fundamental result of this paper that reveals the boundedness of the second 
moments for SDEs (|l.ip under Assumptions 12.1 1 and [ 



Theorem 3.6. Let Assumvtions \2.1[ \3.1\ \3.!^ hold, and 9 > 0.5. Then, for any T > 0, there exists a 
constant C(T) > 0, such that the 9-EM scheme p.ip has the following property 

sup sup E||Af < C(r). 

At<Af 0<tk<T 

Proof. By definition p.2p of function F, we can represent the 0-EM scheme p.ip as 

F{Xt,^,) = F{Xt, ) + /(Ai, ) At + .9(A, J A«;*, . 

Consequently writing {F{Xt^^^), F{Xt^^^)) — ||i^(At^^J||^ and utilizing Assumption 12 . II we obtain 

\\F{X,,^,)r =||F(A,J||2 + ||/(A,JfAt2 + |1.9(A,J|pAt (3.5) 
+ 2(^^(AiJ,/(AiJ)Ai + AAfi,^, 

=\\nxt,)f 

+ (2(A,„/(A,J) + |lg(A,jf) M 
+ (l-2^^)||/(AiJ||2A^2 + AMt,^„ 

where 

AM*,^, = ||g(AiJAu;,,^J|2 - ||.g(A,J||2At + 2(F(A,J,5(A,JAii;4,^J 
+ 2(/(AtJAt,5(A,JAw;t,^J, 

is a local martingale. By Assumption 12.11 and the fact that 9 > 0.5, 

\\F{Xt,^,W < ll^(^tJII' + 2aAt + 2/3||AiJpAt + AM*,^,. (3.6) 
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Let N be any nonnegative integer such that NAt < T. Summing up both sides of inequahty p.6p from 
A: = to iV A Am , we get 



NAX„ 



fe=0 k=Q 
N N 

< \\F{Xt,)\\^ + 2aT + 2pY.\\Xt,^,J'At + Y,^Mt,^^[0A^]ik)- (3.7) 



A:=0 



fc=0 



Applying Lemma l3.3[ Assumption [33] and noting that Xt^ and l[Q x^-^{k) are 7^4^. -measurable while Awt^. 
is independent of J^t^ j we can take expectations on both sides of p.7p to get 



n\F{Xt,,,^^jr < mXtoW + 2aT + 2/3E 



N 



.k=0 



By Lemma 



E||-F(^t„..„ , Jll' <\\F{Xt,)\\^ + {2a + 2/3(1 - 2l39At)-^ 2eaAt){T + At) 



+ 2/3(1 - 2peAt)-^E 



N 



,*:=0 



By the discrete Gronwall inequality 

E||F(Xt„,^^^^+Jf < [\\F{Xt„)\\^ + {2a + 2l3{l - 2l39At)-^ 29aAt){T + At)] exp (2/3 (1 - 2(3dAt)-\T + At)) 

(3.8) 

where we use the fact that NAt < T. Thus, letting m — >■ oo in (|3.8p and applying Fatou's lemma, we 
get 

E\\F{Xt^^,)f < [\\F{Xt,)\\^ + {2a + 2/3(1 - 2p9At)-^ 29aAt){T + At)] exp (2/3 (1 - 2l39At)-\T + At)) . 



□ 



By Lemma 13. 2[ the proof is complete. 



4 Forward-Backward Euler-Maruyama Scheme 



We find in our analysis that it is convenient to work with a continuous extension of a numerical method. 
This continuous extension enables us to use the powerful continuous-time stochastic analysis in order 
to formulate theorems on numerical approximations. We find it particularly useful in the proof of 
forthcoming Theorem 14.21 Let us define 

7j{t):^tk, for t e[tk,tk+i), k>0, 

■q+{t) ifc+i, for t G [tk,tk+i), k>0. 



One possible continuous version of the is given by 



X{t) = Xt, 



/(X,^(,))d5+(1 



f{X^^,))ds+ / g{X^(^,))dw{s), t>0. 
Jo 



(4.1) 
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Unfortunately, this X{t) is not J-j-adapted whence it does not meet the fundamental requirement in the 
classical stochastic analysis. To avoid using Malliavin calculus, we introduce a new numerical scheme, 
which we call the Forward- Backward Euler-Maruyama (FBEM) scheme: Once we compute the discrete 
values Xt^ from the 6'-EM, that is 

Xt, = Xt,_, + efiXtJAt + (1 - 9)f{Xt,_JAt + g{Xt,_JAwt,_,, 

we define the discrete FBEM by 

K^.^K+fiXtJAt + giXtJAwt,, (4.2) 

where Xto = Xtg ~ xq, and the continuous FBEM by 

X{t)^Xt,+ f /(X^(,))ds+ / g(X^(,))du.(s), t>0. (4.3) 
Jo Jo 

Note that the continuous and discrete BFEM schemes coincide at the gridpoints; that is, X{tk) = Xt,. 



4.1 Strong Convergence on the Compact Domain 

It this section we prove the strong convergence theorem. We begin by showing that both schemes of the 
FBEM (|4.2p and the 6'-EM (jS.ip stay close to each other on a compact domain. Then we estimate the 
probability that both continuous FBEM (j4.3|) and 9-EM p.ip will not explode on a finite time interval. 

Lemma 4.1. Let Assumvtions \2. ll[KT[\3.5\ hold. andO > 0.5. Then for any integer p > 2 andm > \\xq\\, 
there exists a constant C{m,p) such that 



E 



\\Xt, - Xtjn[o,x^^ik) <Cim,p)AtP, VfceN, 



and for F(x) = x — Of{x)At we have 

\\Xt,f >\\\F{Xt,W ~\\ef{x,)Atf yk^n. 

Proof. Summing up both schemes of the FBEM (14.21) and the 0-EM (|3.ip . respectively, we obtain 

Xt,-Xt,=6{f{Xt,,)- f{Xt,))At. 

By Holder's inequality. Lemma [3.31 and Assumption 13.51 we then see easily that there exists a constant 
C{m,p) > 0, such that 



E 



||Xt„-Xt„ri[o.A,„](iV) -0E[||/(XtJAi-/(Xt„)Airi[o,A„](A^)] <C{m,p)AtP, (4.4) 



as required. Next, using inequality 2|a||5| < e|ap + e with e = 0.5 we arrive at 

\\XtJ^^\\Xt,-9f{Xt,)At + 9f{Xt,)At\\^ > i\\FiXt,)\\-\\efiXt,)At\\)' 
>l\\FiXt,)r-\\efiXt„)Atr. 



□ 
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The following Theorem provides us with a similar estimate for the distribution of the first passage time 
for the continuous FBEM (|43)) and ^-EM jaj]) as we have obtained for the SDKs (fTTI) in Theorem [2?2] 
We will use this estimate in the proof of forthcoming Theorem 14.41 

Theorem 4.2. Let Assumvtions lKTi \3.1l COl hold, and 9 > 0.5. Then, for any given e > 0, there exists 
a positive integer Nq such that for every m > Nq, we can find a positive number Ato = Aio(m) so that 
whenever At < Ato, 

P(??m <T)<e, for T > 0, 
where 'dm = i^i{t > : \\X{t)\\ > m or > m}. 

Proof. The proof is given in the Appendix. □ 



4.2 Strong Convergence on the Whole Domain 

In this section we present the main theorem of this paper, the strong convergence of the 6-EM p.ip to 
the solution of (II. First, we will show that the continuous FBEM (|4.3|) converges to the true solution 
on the compact domain. This, together with Theorem 14.21 will enable us to extend convergence to the 
whole domain. Let us define the stopping time 



where t„i and dm are defined in Theorems 12.21 and 14.21 respectively. 

Lemma 4.3. Let Assumytions [K7[ \3.1l \3. 5\ hold, and 9 > 0.5. For sufficiently large m, there exists a 
positive constant C{T,m), such that 



E 



sup \\X{tAem)~x{tA9m)\\ 

a<t<T 



< C{T, m)At. 



(4.5) 



Proof. The proof is given in the Appendix. 



□ 



We are now ready to prove the strong convergence of the 9-F,M p.ip to the true solution of 

Theorem 4.4. Let Assumytions \2.1[ \3.1[ U75\ hold, and 9 > 0.5. For any given T = N At > and 
s G [1,2), 9-EM scheme (j3.ip has the property 



lim m\XT~x{TW =0. 

At^O 



(4.6) 



Proof. Let 



Applying Young's inequality 



3(r) = - x{T). 



5s ^ 2 — s 

X y < — X H 

2 2^5^ 



Va;,y,(5 > 0, 
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leads us to 



E||e(T)||^ = E[||e(r)||n{,„>T^^,„>T}] +E[||e(T)||n{,„<T or ^,„<T}] 

]+n\\XT-x{T)ri{r 



+ [||e(T)|p] + ^-^P(t„ < T or i?,„ < T). 

2 2^2-. 

First, let us observe that by Lemma WA\ we obtain 

n\XT - l(r)|in{,„,>T,^„>T}] < Cljn, s)At\ 
Given an e > 0, by Holder's inequality and Theorems 12.21 and 13.61 . we choose S such that 

< <5.E [||x(T)f + II^Tf ] < |. 
Now by (|2.3p there exists A^o such that for m > Nq 

2^2-= 3 

and finally by Lemmas 14. 1[ 14.31 and Theorem 14.21 we choose At sufficiently small, such that 



2' 



s-l 



E[||l(r) - x(T)||n^,„>i.,^,„>^}] + E[\\Xt - XiTWM^^^T.^^^r}]] + ^n^m < T) < |, 



2S^ ^ - ' - 3' 

which completes the proof. □ 



Theorem 14.41 covers many highly non-linear SDEs, though it might be computationally expensive to find 
the inverse of the function F{x) = x — df{x)At. For example, lets consider equation (|1.4p with 
— a + sin(a;)^, a > 0, that is 

dx{t) = (a + sin(a;)2 - ax{tf)dt + (3x{tfdw{t), (4.7) 

where a, /3 > 0. This type of SDE could be used to model electricity prices where we need to account for 
a seasonality pattern, |31) . In this situation, it is useful to split the drift coefficient in two parts, that is 

f{x)^h{x) + h{x). (4.8) 

This allows us to introduce partial implicitness in the numerical scheme. In the case of (|4.7p we would 
take fi{x) = —ax{tY and f2{x) = a + sin(x)^. Then a new partially implicit scheme has the 

following form 

= Yt, + 6h{Yt,^,)At + (1 - e)h{Y,,)At + /2(ytJAi + g{Yt,)Awt,. (4.9) 

It is enough that /i satisfies Assumption 13.11 in order for scheme (14. 9p to be well defined. Its solution 
can be represented as 

Yu^, = {Yt, + (1 - 9)h{Yt,)At + MYtJAt + giYt,)Awt,), 

where 

H{x) = x-efi{x)At, (4.10) 



All results from Sections 3 and 4 hold, once we replace condition (j2.ip in Assumption 12 . 1 1 and Assumption 
[Qbv (|4lI1) (|4l^ ). respectively. 
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Theorem 4.5. Let Assumvtion \3.5\ hold and At G [0, (niax{L, 2/3}0) In addition we assume that 
for x,y € M", there exist constants L,a, /3 > such that 

{x-yJ,{x)-h{y))<L\\x-y\\^ (4.11) 



{x, fix)) + ^\\g{x)r + [(1 - 9){f,{x)J,{x)) + \\\f2{x)r + ^(1 - mihi^W] At<a + /3||x|p. 

(4.12) 

Then for any given T > and s € [li 2), 9 -EM scheme (j4.9p has the following property 

lim E||yT-a;(T)||" =0. (4.13) 

Proof. In order to prove Theorem l4.5l we need to show that results from sections 3 and 4, proved for (13.11) . 
hold for (j4.9l) under modified assumptions. The only significant difference is in the proof of Theorem 
13.61 for (|4.9p . Due to condition (j4.1ip we can show that Lemma [3.21 holds for function H. Then by the 
definition of function H in (|4.10p , we can represent the 6'-EM scheme (|4.9p as 

H{yt.^^) = H{Yt,) + fiYtJAt + g{YtJAwt,. 

Consequently writing {H{Yt^^^), H{Yt^^J) = \\H{Yt^^^)\\'^ we obtain 

\\H{Y,,^,)r = WHiY^l' + WfiYtJW'At^ + ||g(y,,)|pAt (4.14) 
+ 2{H{YtJ,f{YtJ)At + AMt,^, 

= \\H{Y,jr + (2(r*,,/(y,j) + \\g{Ya\') At 

+ [2(1 - 0)(/i(r*j, /2(r.J) + \\f2{Yt,)r + (1 - 20)11/1 (y,jf] At' + am,,^„ 

where 

AMt,^, = \\g{Yt,)Awt,^,f~\\g{YtJfAt + 2{H{YtJ,g{Yt,)Awt,^,) 
+ 2{f{Yt,)At,g{Yt,)Awt,^,}. 

Due to condition (|4.12p we have 

\\HiYt,^jr < \\HiYt,)f + 2aAt + 2/3 ||y,J|2 + AA^,. 

The proof can be completed by analogy to the analysis for the 0-EM scheme (13. ip . Having boundedness 
of moments for (j4.9l) we can show that (j4.13p holds in exactly the same way as for scheme (j3.1l) . □ 



4.3 Numerical Example 



In this section we perform a numerical experiment that confirms our theoretical results. Since Multilevel 
Monte-Carlo simulations provide excellent motivation for our work jl2| . here we consider the measure 
of error (|4.13p with s ~ 2. Although, the case s = 2 is not covered by our analysis, the numerical 
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experiment suggests that Theorem 14.41 still holds. In our numerical experiment, we focus on the error at 
the endpoint T = 1, so we let 

e'll"'^^ ^E\\xiT)-XT\\^. 

We consider the SDE 

dx{t) = ill- ax{tf )dt + I3x^{t)dw{t), 

where {fJ,,a,(3) = (0.5,0.2,70^). The assumptions of Theorem l44l hold. The 9-EM ([3T]) with 6* = 1, 
applied to (|1.4[) writes as 

Xt,^, = X,, + {tJi~ c^XlJAt + pxlAwt, . (4.15) 

Since we employ the BEM to approximate (ll.4p . on each step of the numerical simulation we need to 
find the inverse of the function F{x) = ax^At + x. In this case we can find the inverse function explicitly 
and therefore computational complexity is not increased. Indeed, we observe that it is enough to find 
the real root of the cubic equation 

aXl^^At + Xt,^, - {Xt, + fiAt + pXlAwt,) = 0. 

In Figure[T]we plot e^j°"^ against At on a log-log scale. Error bars representing 95% confidence intervals 
are denoted by circles. Although we do not know the explicit form of the solution to (|1.4p . Theorem 14.41 



Backward EM Scheme 




10"'' 10"' 10"' 

At 



Figure 1: A strong error plot where the dashed line is a reference slope and the continuous line is the 
extrapolation of the error estimates for the BEM scheme. 

guarantees that the BEM (|4.15p strongly converges to the true solution. Therefore, it is reasonable to 
take the BEM with the very small time step At = 2~^^ as a reference solution. We then compare it to 
the BEM evaluated with timesteps (2^ Ai, 2'^ At, 2^ At, 2^At) in order to estimate the rate of convergence. 
Since we are using Monte Carlo method, the sampling error decays with a rate of I/^/M, M- is the 
number of sample paths. We set M = 10000. From Figure[T]we see that there appears to exist a positive 
constant such that 

^strong ^ (J foj- sufficiently small At. 
Hence, our results are consistent with a strong order of convergence of one-half. 



13 



5 Stability Analysis 



In this section we examine the globally almost surely asymptotic stability of the 0-EM scheme p.ip . The 
stability conditions we derive are more related to the mean-square stability, [IHIIID]. First, we give some 
preliminary analysis for the SDEs (|l.ip . We give conditions on the coefficients of the SDEs that 
arc sufficient for the globally almost surely asymptotic stability. Later we prove that the 0-EM scheme 
(I3.ip reproduces this asymptotic behavior very well. 



5.1 Continuous Case 

Here we present a simplified version of the stochastic LaSalle Theorem as proved in ^j, using the 
Lyapunov function V{x) = ||a;|i^. 

Theorem 5.1 (Mao et al.). Let Assumvtion \2.l\ hold. We assume that there exists a Junction z € 
C(E";R+), such that 

{xJ{x)) + \Mx)r<-z{x) (5.1) 
for all X G M". We then have the following assertions: 

• For any x^ G M", the solution x{t) of has the properties that 

limsup ||a;(t)|p < cx) a.s and 
lim z{x{t)) = a.s. 

t—^oo 

What is more, when z{x) — if and only if x — then 

lim x{t) =0 a.s V G M". 



5.2 Almost Sure Stability 

We begin this section with the following Lemma. 

Lemma 5.2. Let Z = {Zn}neN be a nonnegative stochastic process with the Doob decomposition Zn = 
Zq + Al^ — A^ + Mn, where A^ — {Al^}nefi o,'^d A^ = {A^}„gN are a.s. nondecreasing, predictable 
processes with A\ — A^ — 0, and M — {Af„}„£N is local {J>i}„gN-martm(;a/e with AIq — 0. Then 

< lim Al^ < cxD > C lim yl^^ < oo > n lim Z„ exists and is finite > a.s. 

The original lemma can be found in Liptser and Shiryaev [30; . The reader can notice that this lemma 
combines the Doob decomposition and the martingales convergence theorem. Since we use the Lyapunov 
function V{x) = our results extend the mean-square stability for linear systems, Higham [T5J[TB], 

to a highly non-linear setting. The next theorem demonstrates that there exists a discrete counterpart 
of Theorem O for the ^-EM scheme ([3T1) . 

Theorem 5.3. Let Assumptions \2.1i \3.1\ and \3.5\ hold. Assume that there exists a function z G 
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C(IR";R+) such that for all x€W and for all At £ (0, (max{L, 2/3}6')~i), 

ln.._M,2 , (1-20),,.. M,2, 



+ 211.9(^)11 + '-^\\f{x)rM < -z{x). (5.2) 



Then the 9-EM solution defined by iS.l]) . obeys 



limsup ||Xfj^ II < oo a.s., (5.3) 
lim z{XtJ = a.s. (5.4) 



// additionally z{x) — iff x = 0, then 



lim Xt^ = a.s. (5.5) 



Proof. By (|3.5p we have 



||F(X,,,Jf =\\F{Xa\' + (2(X,,, /(X,J) + ||g(X,Jf ) At 

+ (l-2^^)||/(J^iJ||2A^2 + A^v,, 



where 



AMi,^, - ||g(XiJAu;,,^,f - ||5(X,jfAt + 2(F(X,J,5(XtJA«;,,^J 
+ 2(/(XiJAt,.9(X,JAu;,,^,), 

so X^feLo is 9' local martingale due to Assumption 13.51 and Lemma 13.31 Hence, we have obtained 



the decomposition required to apply Lemma 15.21 that is 

||^^(X,„^Jf = \\F{X,,)r -Y^A,,At + J2 AM; 



AT AT 

A Tl ,f 

k- 



k=0 k=0 

where 



A, = - {{2{X,„f{XtJ) + \\9{xa\') + (1 - 2e)\\f{X,jrAt) . (5.6) 



By condition (|5.2p . J2k=o^tk^^ is nondecreasing. Hence by Lemma [Ol we arrive at 

lim ||^^(XiJ||2<oo. (5.7) 



Consequently, by Lemma l3.2l we arrive at 

limsup < oo a.s.. 

By Lemma 

oo oo 

J24XtJAt<Y,At,At<oo a.s, 

/c=0 fc=0 

which implies 

lim z{Xt^) = a.s. 
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This completes the proof of (|5.4p and (|5.5p . 



□ 



A Proof of Lemma 13.31 



Proof. By p.6|) we obtain 



\\FiXt,)f <||F(Xt,_ + 2aAi + 2/3\\Xt,_,f At + AM, 



where 



+ 2(/(X,,_JAt,5(X,,_JAt.iJ. 

Using the basic inequality (ai + 02 + 03 + 04)^*/^ < 4P/^~^(af^^ + Cj^^ + 03^^ + 14^^), where > 0, we 
obtain 



As a consequence 

E [||^^(X,jri[o,A„](fc)] <4f/^-i [\\F{Xt,_J\\n[o,x^]{k)] + (2aAi)P/2 + (2/3m2Atf 



E 



AM,, |P/2l[o,A„](fc) 



Due to Assumption 12 . 1 1 we can bound and ||5(a;)|| for ||x|| < m. Whence, there exists a constant 

C{m,p), such that 



E 



AM,, 



l[o,A„](fc) 



|l5(A,,_JA^,.r + ||g(X,,_jrAt^'/2 + (2||^^(X,,_J|||lg(X,,_JA^,J|)^/2 



(2||/(A,,_JAi||||.g(X,,_JA«;,J|)''/2 



l[0,A,„](fc) 



< C(m,p)E 



l + \\giXitk-i))Awkr 



l[0,A„](fc), 



By Holder's inequality 



E 



AM,, 



l[0,A„](fc) 



< C{m,p) 



1 + (E J f n[o,A„] (A)] ) (E|| Au;„_ Jpf ) 



Hence 



E[||^^(X,jri[o.,„.](fc)] <C(m,p)(l + (E[||.9(A,,_J||n[o,A„](fc)]T^'(lE 
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Since there exists a positive constant C{p), such that E||Awtj._JpP < C{p), we obtain 

E[||F(X,jri[o,A,„](fc)] <Cim,p). 
We conclude the assertion by applying Lemma [321 CH 



B Proof of Theorem 14.2 

Proof. By the Ito formula 



||X(rAi9™)f = llxolP + (2(X(s),/(X^(3))) +trace[/(X„(,))/„x„.9(^^(.))]) ds 



2 / {X{s),g{X^^,)))dwis) 







TAlJ„ 



2 / {X{s),g{X^^s)))dwis) 







TA)?„ 



<||xo|P+ / (2(X„(,),/(X,(,))) + ||g(X,(,))||2)ds 







+ 2 / 

Jo 

+ 2/ (X(s),g(X,(,)))ci«;(s). 
Jo 

By Assumption 12.11 for ||x|| < m 

< 2(||/(x) - /(0)|p + ||/(0)|p) < 2iCim)\\xr + ||/(0)|P), (B.l) 
Mx)r < 2i\\g{x) - 5(0)f + ym') < 2{C{m)\\x\\' + ||g(0)f ), (B.2) 

and 

E\\X{TAArM^<\\xo\\^ + 2aT + 2l3E / - X(s) + X(s)||2ds + C(to)E / - X(s)||ds. 

Jo Jo 

Using the basic inequality (a— 6 + < 2(||a — &|p + ||cp) and the fact that J^^^"^ — X{s)\\'^ds < 

Cirn) Jq^^"^ ll^r;(s) ^ ^{-'')\\d'S, we obtain 

E\\X{T /\-d,nW <\\xo\\^ + 2aT + Al3 [ E\\X{s A i}^)\\^ds 

Jo 

+ 4/3E / -X(s)fds + C(m)E / - l(s)||ds 

Jo Jo 

< \\xof + 2aT + 4/3 / E|lX(s A i?™)f ds + C(to)(4/3 + 1)E / - l(s)|lds. 

Jo Jo 

(B.3) 
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Since \m > a.s., Lemma |4 . 1 1 gives the following bound 



- X^(,)||d.s < C{m,T)M. (B.4) 



To bound the term E " - X{s)\\ds in (|B31) . first we observe that 



Then, by (|RT|) and ((B:2|) 

Il^^(s) - ^(s)llds < C(m,r)At^ (B.5) 



where C(m,T) > is constant. Combining (jB.4l) and (|B.5|) leads us to 

E / - X{s)\\ds < E / - X{s)\\ds 

Jo Jo 

Jo 

< C(m, r)Ati (B.6) 

Therefore 

i-T 

E||X(rAt?„0|P < ||a;o|P + 2ar + C(m,T)At^ +4^ / E||X(s A i?„)|pds. 

"'0 

By Gronwall's inequality 

E\\X{T M)„,)f < [\\xof + 2aT + C{m,T)At^]eyii>{Al3T). (B.7) 

Now we will find the lower bound for ||X(i9„i)|p. From the definition of the stopping time i?„i, if 
inf{t > : ||-^(i)|| > m} < inf{t > : ||-'^^(f)|| > m} then ||X(z9,„)|P = m"^. In the alternative case, 
where inf{t > : > m} > inf{t > : > m}, we have ||X^„Jp > m^. From Lemmas 15^ 

and l4.1l we arrive at 

> i ((1 - 2(39At)\\X^J\^ - 20aAt) - \\ef{xo)At\\\ 
Hence there exist positive constants ci and C2 such that 

||^(^m)||' >Cim2-C2Ai. 

We have 

E|lX(r A 79„0I1' > E [l{^,„<T}ll^(^?m)ll'] > n^,n < T){cW - C2At). 
which implies that 

p., ^ rr., ^ [IkolP + '2aT + Cjm, r)AtV^] exp(4/?r) 

Cim"^ — C2At 

Now, for any given e > 0, we choose A'o such that for any m > Nq 

[llxof + 2ar]exp(4/3r) ^ e 
Cim^ — C2At ~ 2 ' 
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Then, we can choose Aio = ^to{m), such that for any At < Ato 



exp(4/3T)C(m,T)Ati/2 
cim? — C2At 



< 



2' 



whence P(i?m < T) < e as required. 



□ 



C Proof of Lemma 14.3 



Proof. It is useful to observe that since the constant C{T^m) in (14. 5p depends on m, we can prove the 
theorem in a similar way as in the classical case where coefficients / and g in (11.11) obey the global 
Lipschitz condition [THimj. Nevertheless, for completeness of the exposition we present the sketch of 
the proof. For any Ti € [0,r], by Holder's and Burkholder-Davis-Gundy's inequalities 



E 



sup \\x{t^em)-x{t^e^)\\ 

a<t<Ti 



< 2 TE 



- f{x{s))fds + 4E / \\g{X^(s)) - g{x{s))fds 

Ja t 



By Assumption 12.11 there exists a constant C{m) 



E 



sup \\X{tAern)-x{t/\e,n)\\ 
.0<t<Ti 



< 2C(m) TE / \\X^(s)~x{s)fds + 4E \\X,,i^s) - x{s)fds 

\ Jo Jo J 



< 4C(to) TE 



TiAe„ 



\\X{.s)-xi.s)f + \\X^^s)-Xis)\\ 



■4E 



/•Ti A0„ 

Jo 



\\X{s)~x{s)f + \\X^is)~-Xis)\\ 



ds 



Ti 



< 4C(m)(r + 4)E / \\X{s A 0m) - x{s A 0m)\\ ds 
Jo 

rTlABra 

|2, 



+ 4C(m)(r + 4)E 







\\X,,{s)^X{s)rds. 



By the same reasoning which gave estimate (|B.6[) , we can deduce that 



E 



\X^{s) ~ X{s)\Yds < C{m,T)At. 



Hence 



E 



sup \\x{tAe,n)~'x{t/\e,n)\\ 

0<t<Ti 



< 4C(m)(r + 4)C(to, T)At + 4C(m)(r + 4)E 
The statement of the theorem holds by the Gronwall inequality. 



Ti 

sup \\X{t A9m) - x{t Ae^jW^ds 

0<t<s 



□ 
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